This invention relates generally to methods, apparatus, and software for reconstructions of images in computed tomography (CT), and more particularly to methods, apparatus, and software for iterative reconstruction techniques.
Image reconstruction from helical CT scan data benefits from the flexibility of iterative reconstruction techniques. Conventional filtered backprojection suffers from interpolation errors inherent in the approximation of in-plane measurements from helical scans. These problems become especially pronounced with multi-slice detector arrays. Iterative techniques rely on successive approximations of an image, with adjustments of image values based on the difference between measured data and simulated measurements from the candidate image. Let x be a vector containing a 3-D reconstruction of the object; the elements in x will be numbers representing X-ray densities of volumetric elements, called “voxels,” in the three-dimensional object to be imaged. Furthermore, let y be the actual measurements, representing samples from what will subsequently be called the “sinogram,” and let F(x) be the expected values of the sinogram when the 3-D cross section being reconstructed is assumed to be x. Importantly, the model F(x) includes the precise geometry of the helical scan pattern and the source/detector structure, so it can directly account for the helical scan measurements. The difference between the measurements y and their expected values is commonly referred to as “noise,” and may be incorporated into the model in the equation y=F(x)+n, where n represents the noise vector. The optimization problem in iterative reconstruction may be expressed as
                              x          ^                =                  arg          ⁢                                          ⁢                                    min              x                        ⁢                          {                                                                    ∑                                          i                      =                      0                                        M                                    ⁢                                                            D                      i                                        ⁡                                          (                                                                        y                          i                                                ,                                                                              F                            i                                                    ⁡                                                      (                            x                            )                                                                                              )                                                                      +                                  U                  ⁡                                      (                    x                    )                                                              }                                                          (        1        )            F transforms an image x in a manner imitative of the CT system, y is the available sinogram data and the functional Di is a function which penalizes distance between measurement i and the corresponding simulated i-th forward projection of x. U(x) is a regularization term which penalizes local voxel differences. A common embodiment of (1) takes the form
                                          x            ^                    =                      arg            ⁢                                                  ⁢                                          min                x                            ⁢                              {                                                                            ∑                                              i                        =                        0                                            M                                        ⁢                                                                  w                        i                                            ⁢                                                                                                                                                            y                              i                                                        -                                                                                          F                                i                                                            ⁡                                                              (                                x                                )                                                                                                                                                              2                                                                              +                                      U                    ⁡                                          (                      x                      )                                                                      }                                                    ,                            (        2        )            where wi is a constant which weights the contribution of measurement i to the objective function. Frequently, a linear model of the form F(x)=Ax is used, linearizing the relation between x and y with a matrix A. At each iteration, the algorithm finds a perturbation of x which will decrease the value of the above expression. Following adjustment of x, new forward projections by the operator F allow calculation of directions for further improvements. The computation of these simulated forward projections extracts a high computational cost, since storage of the matrix (A) describing the mapping between the image and the sinogram remains infeasible in practical computing systems, and components of F must be recomputed and discarded at each iteration. The quality of the iterative reconstructions depends strongly on the degree to which the mapping mirrors physical reality in the CT scanner.
The basic operation for both forward and backprojection in iterative reconstruction is the computation of the effect of a single element in a digital three-dimensional image on the attenuation measurements of the sinogram. In a typical reconstruction, a given entry in the matrix A of Ax must be used twice as many times as there are iterations, and the number of non-zero entries is extremely large. Because voxels may be modeled as representing rectangular solids in three dimensions, the attenuation of voxel j on sinogram measurement i may be approximated by the length of the line segment between a source and a detector i lying within voxel j, as illustrated in FIG. 3. Such a computational model suffers from several limitations. First, the line-clipping algorithm commonly used to find the length of this line segment consumes the greatest share of iterative computation, saddling these methods with reconstruction times orders of magnitude larger than those of conventional methods. Second, the line segment length is but a coarse approximation of the effect of a volume element on a detector's radiation. The non-zero sizes of both detectors and X-ray source cause blurring of these effects and can limit resolution of reconstructed images. The two-dimensional dosage profile of the X-ray focal spot is often unknown even if its shape is fixed.